{
 "cells": [
  {
   "cell_type": "code",
   "execution_count": 3,
   "metadata": {},
   "outputs": [
    {
     "name": "stdout",
     "output_type": "stream",
     "text": [
      "\n",
      "Created on Sun Apr 19 08:57:13 2015\n",
      "\n",
      "@author: shifeng\n",
      "\n"
     ]
    },
    {
     "data": {
      "image/png": "iVBORw0KGgoAAAANSUhEUgAAAYAAAAEZCAYAAACervI0AAAABHNCSVQICAgIfAhkiAAAAAlwSFlzAAALEgAACxIB0t1+/AAAIABJREFUeJzsnXd4FVX6xz9vIJSQQkIgJEACUqQXxQIKBFCQNYArrBQX0dXfWmFFXFAQAVl1ZVmxr2VdlKq4iC7SlBJaKBZ6UZoQEnoCCTWBvL8/7s313tySm3CTm3I+zzNP7pw55Z3JzLxz5px5v6KqGAwGg6H8EeBvAwwGg8HgH4wDMBgMhnKKcQAGg8FQTjEOwGAwGMopxgEYDAZDOcU4AIPBYCinGAdQihGRwSKyxN92+BsRqSciGSIixdhmnIjkiEiZuIZEZIeIdC5EuTJ7DopIFxFJ9rcdRYmY7wB8g4j8CtQCrgDngKXAk6p6wZ92lUVE5CDwsKqu8KMNccABIFBVc/xlh9WWHKCRqh4o4nbigINARX/vc3EgIl2AGaoa629biooy8fRSQlDgblUNBdoC7YDn/WtS4RCRCuWxbX/hg30u9FNcAdsWa1vF1tMyFC3GAfgWAVDVE1h6AG1tG0QqicgUETkkIkdF5D0RqWy3va+IbBaRsyKyV0R6WNNDReTfIpIqIskiMin3VYeIDBWRNdbf74nIPxyMEflKRJ62/o4Wkf+KyAkR2S8iw+zyjReRL0RkhoicAYY67ZjFjunW8gdFZKzdtqEislZE3haRMyKyS0S65SnraR/WisjrInIKGC8i14nIchE5ZW1vpoiEWvNPB2KBBdbXPs/mfR0jIitF5CVrvRkiskREIuzseUBEfhWRkyLygnV/bPbm2e8qIvJPa/50EVlt938T4I/W/+kJERljV+4mEUmylkmxHpuKdttzROQJEfkF+MWa9oaIHLaeA9+LyO12+QNEZIyI7LPu0/ciUldEVlnt2GZN/4M1f4L1fEq3HodWdnUdFJFRIrIVOCciFeyPgdX27612HBWRKdaiq6x/z1jbusX+HLSWbSEi34rIaWvZ59wc17zXw79yj6vVtg12/8/HRWS7iFSyrs+1lkkXkUQRaW5X7zQReVdEFolIpoisEZEoEZkqImnWc7NNnmPxnIjstNr8cW47Lmx2ew2VWlTVLD5YsHSNu1l/1wW2Aa/bbZ8KfAWEAdWAr4GXrdtuBs7YlY8Gmlh/zwfeA6oAkcAG4P+s24YCq62/OwGH7NqrDlwAorDcIH4AxgIVgPrAPuBOa97xwGWgt3W9sov9m261JQiIA34GHrKzIxsYbq3/Puv+VPdyH7KBJ7A8kFQGGgLdgYpADSAxz7E8CHS1W48DrgIB1vWVwF5rPZWt669YtzUHMoEO1vr/Yd33bm7+r+8CK4Da1uN4KxBobTMH+ACoBLQGLgHXW8vdYP2/ChaHtRMYbldvDpaHhLDc4w0Mtv7fAoARwFGgknXbX4GtWF71ALQCwu3qamBXdzvgONDe2v4Q6zELtDt+PwExdm3bn79JwP3W30HAzXmOs9i1ZX8OBgOpwNPWY1INuMnNcfV0PYj1f/4i0AhIA1rblX3Qalcg8Dqw2W7bNOAEloevSsByLK/q7rfWOwlYkedc2mY9FtWBtcBL1m1dgMN2Nrm9hkrr4ncDyspiPZEyrEsO8B0Qarf9XJ6LtANwwPr7feCfLuqsheWmUtkubWDuCWx/8VnXfwVut/5+BFhm/X0L8Gueup8DPrb+Hg8keti3ACw3yevt0v6cx44jecpstF503uzDr+7atubpC/yY51h3s1t35QDG2G1/HFhk/T0OmGW3rSpuHID1or8AtHSxLbfN6Dz7fJ+bffgLMM9uPQfoks9+pwGtrL/3AAlu8uUA19mtvwdMzJNnD9DJ7vgNdXH+5jqAROs5UcPNPgfYpdk7gIH2/6d89s3t9WDX1mlgFzDKQz3VrfsfYl2fBnxgt/0pYKfdeksgLc9+/5/dei9gr/W3vQPweA2V1sXWJTX4hL6qulJEOgGzsTztZohITSxPLD/KbxNVAvjtXWo9YKGL+uKwPOUctZYT63LYTfufA4OwPMUMBmZY02OBOiKSZl0Xa/ur7cp6mu0QieVp2b7dQ0Adu/WUPGUOYXmq8mYfHNoWkVrAm1h6NcFYnrjSKBjH7H5fsNaD1SZbe6p6UUROu6kjEksPwtPg6nFX7YhIYyxPp+2xOJmKwI95yh6xXxGRZ4E/YekBAoRYbQDLOeLtIG8c8IDdKwrB8j+Icdd2Hh7G8qS8R0QOYHkidnV+5qUesD+/TF5cD6jqIRFZieWG/J5d2QDgFaA/lmOj1iUSS88OHP8nF12sB+OI/bHIPW/z4s01VOowYwC+JXcMYA3wKfBPa/opLDeHFqoaYV2qq2qYdXsyltcVeUnG8vRcw1om3FqutZv25wD9RSQWyxPLPLt6Dti1Ha6qYara266setivU1he08TZpcXheNOvgyOxWF4HeLMPedt+BctTXQtVrQ78EceBR0+25sdRLK/oABCRqlheM7nilNV2V/+b/PgXsBtoaN2HsTgPntr2w/q+/69Af+sxCsfSm8wt4+4ccUUyltcp9v/vYFX93FXbeVHV/ao6WFVrApOB/1qPU37H3Vsb87seEJG7sfQKlgNT7MoOBnpj6a1Ux/IqJvehorDUs/sdh+W8zYs311CpwziAouMN4E4RaaWW/uJHwBvWpx9EpI5YB3qBj4GHRKSrWIgRketV9RjwLTBVREKs264TN/O1VXULlm7zv4Elqpph3bQJyLQOrlWxDvq1EJH23uyIWqb8zQVeFpFgsUwHHMFvPQyAWiIyTEQqWgcim2J57VKgfbASguUVQaaI1MFyY7TnGHBdnjRvbwD/BXqLyK0iEghMcJfR+n/7D/C6dQAwwK5cfm2GABmqekFEmmJ5DeWJECxO9rR1gPRFa1ou/wYmiUgjABFpJSLh1m15j8dHwGMicrM1bzUR+Z2IVMvHBqz57xeR3J7HWSw3/hzgpPWvu5v8N0BtERlu3YfgXBvsye96sLb9EZbe0INAgoj0shYPwfLKLt26P69S8AeCvP+3J63tRwBjgM9clLmma6ikYhyA73A4CVX1FJZewIvWpOewDBptEMtMm2+BJta83wMPYXEaZ7G8g82de/wAlsGsXVheg3yBZUDSHbOxDKDOsrMlB0jAMjB2EMsg2UdAaAH2bziWp7YDWLq9M1V1mt32jUBjLE93k4B+qppeyH2YCNyIZSB5Ab/1ZHL5OzDOOqvjmdzdtNvu6el2FzAMy+uyVCxP2Sew3FRc8SywHfgei3P9O79dN3nb0Tzl7heRDCwDxXlvKnnLLrUuv2D5H13A8dXY61ic8LcichaLQ6hq3TYRmG49Hv1V9Ufg/4B3rK8sfsFxZper42Ofdhew02r7VGCAql5W1YvAy8A6a1sON3dVPQfcCfTB4pR+AeJdtAUwGjfXA5bjNV9Vl6pqGpbxrI+sDm86lteHKcAOLAPWBSXv/s+2tr8Py+SBl50K+OYaKnGYD8EM14yIDMXyYVaBvyT1N9anyDNYZtcc8rc9huJFSsBHhf7E9AAM5Q6xzJGvar35/xPYZm7+hvKIcQCG8khfLK9/jmB5nz3Qv+YY/Ei5fgViXgEZDAZDOcX0AAwGg6GcUmo+BBMR01UxGAyGQqCqLqcsl6oeQGE/dx4/frzfP7ku7sXsc/lYzD6Xj+Va9tkTpcoBGAwGg8F3GAdgMBgM5ZRy4QDi4+P9bUKxY/a5fGD2uXxQVPtcpNNAReRjLJ9PH1c3AcxE5C0sEf/OAw+qJZ6Nq3xalLYaDAZDWUREUD8NAk8DerrbaA3w1FBVGwOPYomLbzAYDIZioEgdgKquBdI9ZOmLJbgTqroRCBORqKK0yWAwGAwW/P0dQB0cIx6mWNOOu85uMBgMxcORI0fYu3evU3pwcDA33XSTU3pmZiY//PBDkeUvCvztAAwGg6FE8sknnzBu3Din9DZt2rBli/NQ5f79++nWrZtP82dkZHD8+HEaN25cmF3IF387gBQc1Xjq4iwtaGPChAm23/Hx8aV+NsDatRFcueLpDZnBYPAXHTvCoEGwe7djekzMVhITncdUU1KgbVvnegqTf+XKAM6c6cLp072pUWMhKSkriI/3bhJMYmIiiYmJXuUt8mBwIlIfWKCqrVxs+x3wpKreLSK3Am+o6q1u6ilzs4ASE8Xrf6qh6BCBazq1rrmC8ockJqIl8AEuURKJ13jbukwUdHzx/m8zMjJYtWoVOTk5dOnSherVq19TfZ5mARVpD0BEZmNRBKohIoeB8ViUoVRVP1TVRVapun1YpoE+VJT2GAwGQ0nm5MmTLFq0iHbt2tGyZUsCAop2omaROgBVHexFnqeK0gaDwWDwhqysLH+bQI0aNbjnnnsICwsrlvbKxZfABoPB4IkNGzbQpEkTNrPZr3YEBAQU280fjAMwGAzlnLNnzzJo0CAOHTrEBjYUW7s5OTnF1pY7/D0LyGAwlFOOHDnC7t27CQ8Pp3379k7b09LS+PHHH53SfZ3//fff59dff+XGG2/kkR8fKeTeeI+qsmvXLnbt2kW/fv2K/D2/J4wDMBgMxc+uXTTs2ZOsrCzi4+NZuXKlU5atW7fSo0cPp/SiyF+tWjXmzJlDShO3s9B9gv0MnzvvvNOvN38wDsBgMBQzZ8+ehUmTyMrKom3bttx4440u84WHh3PHHXc4pbdp08an+StVqsQzzzxD48aNSXH/GdI1kfvU/+OPP9K2bdtimeHjDaVGFN58B2AoKsx3AMWHqjJo0CA+//xzbrzxRpKSkqhUqZK/zbJRVN8BnDp1inXr1vlkXn9B8dt3AAaDwf/k5ORw8eJFp3QRISgoqFjz79+/nwULFkCVKsyZM6dE3fyLksjISPr06YOIy/uw3zAOwGAo4xw7dow6deo4pUdHR5Oamlqs+Rs1asSPP/5Isy+/LLL4NiWVknbzB+MADIYyj7sn96pVq/olf9OmTS2Bdsogqsrx48epXbu2v03xCuMADIYyTnR0NOfPny8x+csq9jN8evfuXSIGefPDOACDwWC4BkrqDB9vMA7AYChjLFy4EIAePXoQGBjoZ2vKNufOnWPlypXk5OTQp0+fYp/hc62UDjdlMBi8ZsyYMSQkJNgcgaHoqFixIvXr16d3796l7uYPpgdgMJQptm3bxrZt2wgPD6dXr17+NqfMU6VKFVq1cpI6KTWYHoDBUIaYMWMGAAMGDKBy5cp+tsZQ0jE9AIPBB3z22WesXr0agJiYGF544QWnPCkpKbz88stO6b7Kf/jwYaZMmQLAkCFDCrwPBvdkZGSwfft2OnToUGoGeL3BOACD4Ro5Ctx///228L6tWrVyeYM+ffo0//rXv5zSfZU/Pd2iL92kSRM6dOhQ0N0wuEAQdu7caZvhU9YwDsBguEZmYwmH0LFjR+6//35q1KjhMl9MTAzvvvuuU7qv8tepU4f33nuPO++8s0R+dVrayMjIYETECPbt21cqZ/h4gwkG50dMMLiSgTwXAVXTC12+9QTYBjAAaOYjo8o6XVbCqq7+tsKJlRNW0nVCV2IqxvBMxDOsvrya+cPnl+rXPp6CwRkH4EeMAygZXEvER1VlZUAAcx99lDfffNMMvHqJJCai8fH+NsOJ3GigqkpmZiahoaH+NumaMdFADYYiQkToBnR7/31/m2LwISJSJm7++WF6AH7E9AAKxtqItVxJv+JvMwxljJyKOQRcsbziqRhekdvTbvezRb7F9AAMZYIr6VccxDp8xTWLfpRgQZiS+qqlJJAbw2f79u3079+fihXL3+2w/O2xwWDHiRMnIA0OHDhAlSpViImJccpz8eJFjh496pTuLr+h5GMfufOuu+4qlzd/MA7AUI5ZsGABffr0AaDhWw3p2rUrK1ascMq3fv16unfv7pTuLr+h5FKaI3cWBcYBGMotn332meVHNWhQqwHR0dEu81WpUoUGDRo4pbvLbyi5nD17lgMHDpTZef0FxTgAQ7kkJyeH7777zrLyIBx454DbvB07duTAAffbDaWH6tWrk5CQYD6Us1J++z6Gcs3ly5d54okn6Nu3L0T62xpDcWJu/r9hegCGcknVqlWZMGECYJkFZChbqCopKSnUrVvX36aUaIwDMBgMZQr7GT61a9cutzN8vMEcGYPBUCYwM3wKjnEABoOh1HP+/HlWrFhRarV5/UWRu0cRuUtE9ojILyIy2sX2eiKyQkR+EpEtImJ07AwGQ4EIDAzkuuuuK7XavP6iSB2AiAQA7wA9gRbAIBFpmifbC8DnqnoDMAh4ryhtMhj69evHyJEjOXPmjL9NMfiISpUq0aJFC/PKp4AU9dG6GdirqodUNRv4DOibJ08OkBt2rzqQUsQ2GcoxKSkpfPnll3zwwQcEBQX52xyDwa8UtQOoAyTbrR+xptkzERgiIsnAN8CwIrbJUAYYMWIEjRo1clr+97//ecx/0003ARAfH0+lSpWK02SDD8id4XPliokK6wtKwiDwIGCaqk4VkVuBmVheFzmRO28bLBdwvIlyWG45fvw4+/fvd0o/d+6cV/kHDRpUZLYZfI+q2rR527VrZ171eCAxMZHExESv8hapHoD1hj5BVe+yrj8HqKq+ZpdnB9BTVVOs6/uBW1T1VJ66jB5AOSdXrQng2LFjLm/2UVFRhISEOKXb5w8KCnKI4mnCQZds7Of1d+nSxQzyFhB/6gF8DzQSkTjgKDAQyxO/PYeAO4BPRaQZUDnvzd9gyEvt2rWLNL+hZHD27Fm++uor2rVrZ+b1FwFF6gBU9aqIPAV8i2W84WNV3S0iE4HvVfUb4FngIxEZgWVAeGhR2mQwGEoPoaGh9OvXj+DgYH+bUiYp8jEAVV0CXJ8nbbzd791A2dJgMxgMPkFEzM2/CCkJg8CGckrA8xFolXSv889mNjJAIA6o5js75FK47yozFJorV66YuD3FjDnaBr+hVdILNPj6+ITHYS488sgjfDT5oyK0zFCc5M7w2bp1K/fddx+BgYH+NqncYByAodSwjW0AdO7c2c+WGHyF/Qyfu+++29z8ixnjAAylgpycHOMAyhB55/WbGT7+wTgAQ6lgx44dnOMcsbGxxMXF+dscwzVy7tw5Dh06RN++fc28fj9iHIChVLB69WoAOh8+bPnwqiQRbgaRC0pISAh33323v80o9xgHYPArhw4d4uzZs07psbGxDk+GrVq1ohe96PvFn6B/f5d11a9fn0OHDhWZrW5JTy95TsmOkmuZwZfExcXx66+/FqiMcQAGvzJy5EjmzZvnlP7FF1/Q3+5G36VLFxQlvn+827oOHTpEWQsXYjB4S2HE7o0DMPiV2NhYWrVq5ZQeFhbmB2sMhvKFcQAGv/L666/72wSDodxi5l0ZDAZDOcU4AIPBYCineOUARKSSiDQqamMMZZ/z588TExNDr169wIzXFjlJSUk0adKE0NBQt2ppuUycOJEhQ4a43d6gQQNWrFjhE7uysrJo0aIFx48f90l9ZZ1bbrmF3bt3+7zefB2AiNwNbAe+s663FZH5PrfEUC74+eefOXr0KMnJyeVqfmL9+vUJCgoiNDSUmJgYHnroIS5cuOCQJykpie7duxMaGkp4eDh9+/Z1uugzMzN5+umniYuLIzQ0lMaNG/PMM8+Qlpbmst0XX3yR4cOHk5GRQZ8+ffK1szAzSXIZPXo0kZGR1KxZk+eee85j3g8//JAuXboQFRVV6PZKAlOnTiU6Oprq1avzyCOPkJ2d7Tbv3Llzad68OWFhYbRs2ZKvv/7atm3nzp3cdddd1KxZkwoVKjiV/etf/8q4ceN8br83PYCXgFuAMwCqugUwvQFDodi1axcAzZs397MlxYuIsHDhQjIyMtiyZQubN2/m1VdftW1fv349PXv25Pe//z1Hjx7l4MGDtG7dmttuu802tzs7O5tu3bqxe/duvv32WzIyMli/fj01atRg06ZNLts9dOhQsRzrDz74gP/9739s376dbdu2sWDBAj788EO3+d9//32PvQ1PXL16tbBm+pSlS5cyefJkVq5cyaFDh9i/fz/jx493mTc1NZUhQ4bwxhtvcPbsWSZPnszgwYM5dcqifRUYGMiAAQP4z3/+47J87969WblyJSdOnPDtTqiqxwXYYP272S5tW37lfL1YTC1brFxZ9vYpP55//nkFdPz48cqEgu3/SlZ63F6Sz5H69evr8uXLbeujRo3ShIQE23qnTp30qaeecirXq1cvHTp0qKqqfvTRR1q7dm29cOGCV202bNhQK1SooFWrVtWQkBDNysrS1NRU7dOnj0ZERGjjxo31o48+suWfMGGCDhkyxLY+ffp0jYuL08jISH355Zed9sGejh07OtT18ccfa4cOHVzmPXz4sAYFBenVq1dtaQsXLtR27dppaGioxsbG6oQJE2zbfv31VxUR/fjjjzU2Nla7dOmiqqrr16/Xjh07avXq1bVt27aamJhoKzNt2jRt1qyZhoSEaMOGDfWDDz7w6pgVhMGDB+vYsWNt68uXL9fatWu7zLtx40aNiopySKtZs6Zu2LDBIW3fvn0aEBDgso4777xTp0+f7tYed+e/Nd3lfdWbHsBuEbkPCBCRBiIyFdjgWzdkKC/k9gCaNWvmZ0v8x5EjR1i8eDGNGzcG4OLFiyQlJTl8+JbLfffdx3fffQfA8uXLueuuu6hatapX7ezbt4969erZeh65T5mxsbEcO3aML774gjFjxrgUEN+1axdPPPEEs2bNIjU1ldOnT5OSkuK2rZ07d9KmTRvbeps2bdi5c6fLvNu3b+e6665zCP4WHBzMjBkzOHv2LAsXLuT99993GrNYvXo1e/bsYenSpaSmppKQkMCLL75Ieno6U6ZMoV+/fpw+fRqwaEMvWrSIjIwMpk2bxogRI9iyZYtLe9atW0d4eDgRERGEh4c7/I6IiCApKcnrfT5x4gTp6c4aF+3bt6dZs2YsWLCAnJwcvvrqK6pUqULr1q3dHFFnmjVrxtatW73O7w3eOICngBuxyDV+CVwG/uJTKwzlhj179gD+eQUk4pulsNxzzz2EhoYSGxtLVFQUEyZMACAtLY2cnByio6OdykRHR9teE5w+fdplnvxQ69fRR44cYf369bz22msEBgbSpk0bHnnkEaZPn+5UZt68efTu3ZvbbruNwMBAJk2a5HF84Ny5cw4f74WFhXHu3DmXec+cOUNISIhDWufOnWnRogUALVu2ZODAgaxatcq2XUSYOHEiVatWpXLlysycOZO7776bnj17AtC9e3fat2/PokWLAOjVqxf169cHoFOnTvTo0YM1a9a4tOe2224jPT2dtLQ00tPTHX6npaXRsWNHr/dZVcnMzHTKGxAQwJAhQxg8eDCVK1fmj3/8Ix988IHXzhws8ZPOnDnjdX5v8MYB9FTV0arazro8B/TyqRWGcsOWLVvYtm0bTZs2Lfa2VX2zFJavv/7aFv9+z549tht7eHg4AQEBHD161KnM0aNHiYyMBKBGjRou83hLamoqERERBAUF2dLi4uJcPtmnpqZSr14923pQUBA1atRwW3dwcDAZGRm29YyMDLdSjuHh4U43yU2bNtGtWzdq1apF9erV+eCDD2zHJ5e6devafh86dIi5c+cSERFhe1pft26d7fgsXryYDh06UKNGDcLDw1m8eLFTfdeKq30WESfnBrBs2TJGjRrF6tWryc7OJjExkYcffpht27Z53V5mZqbPI6d64wBecJE21qdWGMoNVapUoVWrVuVS+CP3SbxTp04MHTqUkSNHApaba4cOHfjiiy+cysydO5c77rgDgDvuuIOlS5dy8eLFQrUfExNDWloa58+ft6UdPnyYOnXqOOWNjo62zNSycuHCBdvrFVe0aNHC4fXEli1bbE/0eWndujUHDx4kJyfHljZ48GDuueceUlJSOHPmDI8++qhTXCf7Hki9evV44IEHSEtLsz2tZ2ZmMmrUKLKysujfvz+jRo3i5MmTpKen06tXL7dxotauXUtISAihoaEOS27aunXrvN7nqKgowl1Eh926dStdunShXbt2gOWV0C233MKyZctc1u2K3bt3O7xy8gnuBgeAnsBU4Djwut3yb+B7d+WKaqEED/AVlvI4CGxPeR4EPnnypFarVk23bdumqqpr167V4OBgffvttzUzM1PT0tJ07NixGh4ervv27VNV1cuXL+vNN9+svXr10j179mhOTo6eOnVKX3nlFV28eLFX7Xbu3FmHDRumly5d0q1bt2pUVJSuWLFCVR0HgXfu3KkhISG6bt06zcrK0pEjR2pgYKDbQeD3339fmzdvrikpKZqSkqItWrTQDz/80O3xaNOmja5fv962HhUVZRvg3Lhxo9aqVctmS+4gsP2gcXJyskZHR+vSpUv16tWrevHiRU1MTNSUlBTNzMzUihUr6urVq1VVddGiRRoUFKTjxo1za09hWLJkiUZHR+uuXbs0PT1du3XrpmPGjHGZd9WqVVqzZk3dsmWLqqr+9NNPGhkZqd99950tz6VLl3Tnzp0qInrp0iW9fPmyw7aIiAg9evSoW3vcnf94GAT2dMNtBzwMHLL+zV3uAyLdlSuqpSRf3IXFOIDy4wAaNGjgdPN84okntH///rb1devWaXx8vAYHB2tYWJgmJCTorl27HMpkZGToiBEjtF69ehoSEqKNGjXSkSNHalpamlftpqSkaEJCgkZERGijRo0cbtKuZgHFxsZqZGSkvvLKKy73wZ7Ro0drRESE1qhRQ5977jmPx+O9997Txx9/3LY+b948jYuL09DQUO3du7cOGzbMwQEEBAQ4OABV1U2bNmmXLl00IiJCa9WqpQkJCZqcnGyrPyoqSsPDw/WBBx7QQYMG+dwBqKpOnTpVo6KiNCwsTB9++GHNysqybWvRooXOnj3btv7uu+9qo0aNNDQ0VBs2bKhTp061bct1cgEBARoQEKAiog0aNLBtnzt3rvbr18+jLYVxAKL5vNQUkSqqesmHnY5CISKan62ljcREIT6+bO1TQZCJUiBR+ERJJF7j3dcn4rabbyhZZGVlccMNN7B8+fJS/zHgxJcjAAAgAElEQVRYcdChQwc+/vhjj5Mn3J3/1nSXI/jeRAOtIyIvA82BKrmJqtrEi7IGg42zZ8+aMM8GACpVqsSOHTv8bUapYf369UVSrzeDwJ8A07B8uN8LmAt8XiTWGMos2dnZREZGEh0dTVZWlr/NMRgMeNcDCFLVpSIyRVX3Ay+IyA+A7wNTGEoEqamp5OTkOEy7y2XJkiUcOXLEKb1nz54O0wbz5j916hRXrlyhatWqVKpUqUjsNhgMBcMbB3BZRAKA/SLyGJACOE90NZQJjh49Stu2bfnkk09cOoC33nqLxYsXO6UvXLjQpQPIm79ly5a+NdhgMBQabxzACKAaMBx4GQgD/lSURhn8Q05ODg888AAnT57kp59+4ne/+51Tnp49exITE+OU7spZ5M0fGBjI448/7lujDQZDocnXAajqRuvPTGAIgIg4fzliKPVMmTKFZcuWERkZycMPP+wyz1/+UrAoIAXNbzAYig+PDkBEbgLqAGtV9ZSItABGA90A1498hgKRk5PD3r17bSFuo6KiXH5yf+zYMZcx332Vf9GiRYwZMwaATz/9tFAxZwwGQynD3QcCwKvAL8AXwA/ABOAEMBLLwLD5EOwaWbkSfeyxxxSLNpYC+vrrr7vM+/TTTzvk83X+P//5zwro008/7bP9y4/y9CGYwVDUuDv/KcyHYCKyC7hRVS+KSASQDLRS1QM+9kFecS0fgq2NWMuV9Cs+tsgHrOzK412bsoc91KY2lanMQAZyF3c5ZZ3NbL7lW6d0X+ZPJpkRjKASJXOWTkYI9PWkati1q/kQLA9JSUk8+OCDHDt2jJkzZ3pUBZs4cSL79u1jxowZLrc3aNCAjz/+mG7duvnEtttvv513333X9/FtyhhZWVm0adOGtWvXegzIV5gPwTw9cf+UZ32zu7yeFuAuYA+W3sRoN3nuA3ZikZ6c6SZPoT1jfk+O/mLlSrRJkyYKOH3y71eK6Um6oD2AfOsrwT2AuLg4mzBLdHS0Pvjgg3r+/HmHPOvWrdNu3bppSEiIVq9eXfv06eMyFMRf/vIXjY2NtYWCGDFihJ4+fdplu927d9e3337bKxvzhoLIiydBmJUrV2rXrl01LCzMIYSBOxYsWKC9evXyyq6SzKxZszQuLk6Dg4P197//vaanp7vMt2bNGg0ODtaQkBANCQnR4OBgFRH98ssvnfJ269bNKfbRP/7xDx05cqRHW9yd/xRSEOY6EfnSuswHGtitf+mhnL3nCQDewRJYrgUwSESa5snTCMu4QgdVbQU87U3dZYW7776be++9l5o1a/rbFEMRUtYlIatVq8bDDz/MlClTvMpfFiQhd+7cyWOPPcasWbM4fvw4VatWdTvL7fbbbyczM5OMjAwyMjL45ptvCAkJ4a67HHvjs2fP5sqVK07aC4MGDeLTTz/1qDlcKNx5BqC7p8VduTx13Aostlt/jjy9AOA14E9e1OXR+3miJPcASiSmB+BzyrokZC7Lli3LtweQlZWlVatW1ZSUFFvapk2btEOHDlq9enWNiYnRp556SrOzs23bRUTfffddbdy4sV533XWqqrp792698847NSIiQps2bapz58615fckMekrxowZo/fff79tff/+/VqpUiU9d+5cvmUffPBB/dOf/uSQdvbsWb3++ut148aNLoPfNWnSxBbh1BXuzn8K0wNQ1eWeFi/9Sx0sYwe5HLGm2dMEuF5E1opIkoj09LJug6FUUhYlIQvC3r17qVChgsP3JBUqVOCNN94gLS2N9evXs2LFCt577z2Hcl9//TWbNm1i165dXLhwgR49evDHP/6RU6dO8dlnn/Hkk0/aFOe8kZjMJTk52aMk5GeffeayXF5JyOuuu45KlSrxyy+/eNz/CxcuMG/ePB588EGH9DFjxvDEE0+4DY7XtGlTn0tCevMhWFFTEWgEdAZigdUi0lJVM/JmzJXQA4iPjyc+Pr6YTDSUBWTiNeg52lGQCKb23HPPPYBFSrB79+6FkoRs3759we1VR0nIJUuWOElC5r2W7CUhASZNmsQ777xT4LZd4UoS8oYbbrD9jo2N5c9//jOrVq1i+PDhtvQxY8bYFLHmzp1LgwYNeOCBBwCLHu+9997LF198wbhx4+jcubOtnL3EpKtB8Hr16rnU8c2PvJKQYJGFdCUJac+8efOoWbMmnTp1sqX98MMPJCUl8fbbb3P48GGX5byVhExMTHTp1F1R1A4gBctNPZe61jR7jgAbVDUH+FVEfgEaAz/mrczeARgMBaWwN25f8fXXX9O1a1fWrFnD4MGDOXXqFKGhoQ6SkE2aOAbZLQ5JyB9/dLrUCiwJWRBcSULu3buXZ555hh9++IGLFy9y5coVbrzxRoc8eSUhN2zYQEREBGBxclevXrU5hI0bN/L888+zY8cOsrKyyMrK4g9/+INP7M8lryQkWGQhXUlC2jN9+nSbnbm2P/nkk7z55pseQ5p7KwmZ9+F44sSJbvN6Ew0UABGp7G1eO74HGolInIhUAgYCefthXwFdrW1EYrn5+2WqqcFQlORe2GVRErIgNGrUCFV1cGaPP/44zZo1Y//+/Zw5c4aXX37Z6UaYVxIyPj7eQRIyIyPD1ku5//7785WYzCU5OdmjJOScOXNclssrCXngwAGysrKcnLg9R44cITEx0cEBZGRk8OOPPzJgwACio6O5+eabUVXq1q3rIEdZFJKQ+ToAEblZRLYDe63rbUTkbW8qV9WrwFPAt1imeX6mqrtFZKKIJFjzLAVOi8hOYDnwrKoWvD9WCjl2DKZNm8batWv9bYqhmHn66af57rvv2L59OwB///vf+fTTT3nnnXc4d+4c6enpvPDCC2zYsIEXX3wRgCFDhlCvXj369evHzz//jKpy+vRpXn31VZYsWZJvm3Xr1qVjx448//zzXL58mW3btvHxxx+7nI3Tv39/vvnmG5KSksjOzubFF1/0+I2FqnL58mWysrLIycnh8uXLbmesBAYGcscdd7Bq1SpbWmZmJqGhoQQFBbFnzx7+9a9/edyXhIQEfvnlF2bOnMmVK1fIzs7mhx9+4OeffwYsr2fCw8MJDAxk06ZNzJ49221d9erVc5ihk7vkpg0aNMhlufvvv58FCxawbt06zp8/z/jx4+nXrx/VqlVz29b06dO57bbbaNCggS0tLCyM1NRUtmzZwtatW1m0aBEAP/30E7fccgtg6ZGlp6dz6623ejwuBcbd6HDuAmwA4rD7DgDYkV85Xy+UwVlAEyZYvs699957/W2KI2YWkM8p65KQiYmJDpKGAQEB2rVrV7fHY+HChQ7fAaxevVqbNm2qISEh2rlzZx0/frx26tTJtj0gIED379/vUMcvv/yid999t9asWVMjIyO1e/fuunXrVlX1LDHpS+bMmaOxsbEuvwPo1auXvvrqqw75mzVrptOmTfNYpysJzMmTJxfJdwDe3Hg3Wf/aO4Ct+ZXz9VIWHcCzz1ocwEMPPeRvUxwxDsBQDNx+++02kXSDey5fvqzNmjXTkydPesxXGAfgzSBwsojcDKiIVACGYfmq13CN5L6K9WZgx2Aoa6xZs8bfJpQKKlWqxK5du4qkbm8GgR8HnsEym+c4lo+7TFB3H3DunOWv0ck1GAz+wJsewBVVHVjklpRDch2A6QEYDAZ/4E0P4HsRWSQiQ0XESEFeAxERIPLb0rQpwBCefrqlQ7q/Fyi+dgwGg//I1wGoakPgb8CNwHYR+UpETI+gEKSng2WE1bL06AGq01Ht7pDu7wWKrx2DweA/vPoQTFWTVHU4cAOQAcwqUqsMBoPBUOR48yFYsIjcLyILgE3ASaBjkVtmMBgMhiLFm0HgHcACYLKqmnlbBoPBUEbw5hXQdao6zNz8DYbSR1JSEk2aNCE0NNRtOORcJk6c6FGkpUGDBqxYscJntt1+++0+D29cFsnKyqJZs2Y+i8Vkj1sHICL/tP6cZ68EVhBFMIN7cnJymD8fZs0ywynlgfr16xMUFERoaCgxMTE89NBDXLhwwSFPUlIS3bt3t0UI7du3L7t373bIk5mZydNPP01cXByhoaE0btyYZ555hrS0NJftvvjiiwwfPpyMjAyPesC55FWi8pYpU6bQqlUrQkNDadiwYb7KYN988w2hoaGlXg949uzZ1K9fn5CQEO6991634ZrXrl3rEHAuJCSEgIAA5s+fD8Cnn35KxYoVHYLQrV69GrB8CPbwww87KMj5Ck89gM+tf98B3nWxGK6Bc+fO8dZb8Nhjj/nbFEMxUNYlIQFmzJjBmTNnWLx4Me+88w5z5851m9dIQjpLQnbs2NEhCJ29pkGxS0LmLsBT3qQV9UIZiAVkvwuHDh1SQOvWres/g9xhYgH5nPIiCZnL8OHDdfjw4S63GUlIZ0nITz75xCH4nSuKVRLSjj+5SHvYN+6n/HL27FnAhIEoj5QHScg1a9bQokULl9uMJKRrScjNmzdTq1YtmjZtyt/+9jdycnIctheFJKSnMYABIjIfaJDn/f93QP66ZAaP5L4rNGEgihE/f8J8zz33EBoaSmxsLFFRUYWShHSVJz9UHSUhX3vtNSdJyLzYS0IGBgYyadIkr8cHxo8fj6ry0EMPudzuThLy5ptvRkQcJCHtyZWErFy5Mt98841NElJEHCQhATp37mxzQPaSkK7IlYTMFZax/52WlsbAga6/e/WlJGSXLl3YsWMHJ06cYN68ecyZM4d//OMfDuW8lYQsCJ6mgW4CTmORcbR/558JbPapFcWAeKmRWaSsBEm0/k5KAiBs3bqSFxchPNzfFhQNfv78uDxIQr7zzjvMnDmTtWvXEhgY6DKPkYR0lIQEyySBXFq0aMGLL77IlClTGD16tC3dW0nIguDWAajqQeAgsMynLfoJLQEC8iK/3YM2h4UxfsNY7rzzLRg2zL+GGYqF3Cdxe0nI+fPnO0hCdunSxaFMXknIcePGcfHiRa9fA9ljLwmZq1rlSRIy93UKeCcJ+Z///IfJkyezZs0ajz0Ve0nI3HyPP/44N9xwA59//jlBQUG8+eabzJs3z6GcK0nIpUuXumzj/vvvZ/jw4SxdupTAwEBGjBjh1v7k5GSaN2/u1MNRVUSEDz74wKUq2LVIQn744Ydu89i3b8/u3bt59tln8y1XINwNDgCrrH/TgTS7JR1Ic1euqBbK2CCwqurKlSV30LI4KM+DwCdPntRq1arptm3bVFV17dq1GhwcrG+//bZmZmZqWlqajh07VsPDw3Xfvn2qahEGufnmm7VXr166Z88ezcnJ0VOnTukrr7yiixcv9qrdzp0767Bhw/TSpUu6detWjYqK0hUrVqiq4yDwzp07NSQkRNetW6dZWVk6cuRIDQwMdDsIPHPmTK1du7bu2bPHq+PRt29fnTNnjm395ptv1kmTJqmqZXD3+uuvdxgUFREHRbDMzEytX7++zpgxQ7OzszUrK0u///57W/tRUVE6ffp0VVXduHGj1qpVy+eKYDt37tSwsDBdu3atnjt3Tv/4xz/q4MGDPZZ5+eWXtUuXLk7pixcv1uPHj6uqZf9btmxpOx6qFiW3yMhIzcrKclu3u/OfQg4Cd7X+jQRq2i256waDwUvyPl1GRkYydOhQXnrpJQBuu+02li5dyrx584iOjqZBgwZs3bqVdevW0bBhQ8AyH3zZsmU0bdqUO++8k7CwMG699VZOnz5t047Nr905c+Zw8OBBYmJi6NevH5MmTaJr165O5Zo3b867777LoEGDiImJoUaNGg6vYPIybtw40tLSuOmmm2zz2J944gm3+f/85z87jD1MmTKFWbNmERoayqOPPur03j3vfgQHB/Ptt9/y2WefERMTQ0xMDM899xyXL18G4L333mPcuHGEhYXxt7/9jQEDBri1pbA0b96c999/n8GDB1O7dm3Onz/Pu+/+9rb8d7/7HX//+98dysycOdNp8BcsA/ytW7cmJCSEhIQE+vfvz/PPP2/bPmvWLIYOHer2tVphEc3nvaiI1AdSVTVLRG4HWgMzVTXDY0EfIyKan63uSJRE4jXetwYVAvtXQACJiUJ8fPkNiykTBR3vu/0XEadus6Hk0qlTJ955551S/zFYUZOVlUXbtm1ZvXq1bTzIFe7Of2u6y4FGb2IBfQXcJCINgWnAN8BsIMEb4w0Gg8EVRhLSO/wtCZmjqtnAvcDbqjoCcB41MhgMBkOpwhsHcEVE/gAMwfL0D+DbF1HlkP/+97988w0cO3bM36YYDIZyirdfAnfFEg76gIg0AOYUrVlln3/84x/8859w8OBBf5tiMBjKKfmOAajqDhEZDjQSkabAPlV9uehNK9vkhoIwXwIbDAZ/ka8DEJFOwAwgBRCgtogMUdV1RW1cWSb3k24TC8hgMPgLb2YBTQV+p6q7AESkGRaH0L4oDfM1MrEEhFuYADLR+tv6UWKdd+tAJX8Z5F/Cq5TRkBMGQynBGwdQKffmD6Cqu0Wk1N2yfDnf/Fq5fPkyVSZUoUIFyP5bTqFFOAwGg+Fa8GYQ+CcReV9Ebrcu/6IUBoMrSVy9epWnn36aPn0Kr8BkMHiDkYQsG9xyyy1O6nC+wBsH8BhwABhlXQ4Aj/rcknJEUFAQU6dOZfhwf1tiKC7KuiTkG2+8QcOGDQkLC6Nu3bqMHDnSKZ69PUYS0lES8vPPP6dp06ZUr16d2rVr89BDD3Hu3Dlb+b/+9a+MGzfO5/Z7dAAi0gq4C5ivqn2syz9U9ZLPLTEYyjBlXRKyb9++bN68mbNnz7Jjxw62bNnCW2+95Ta/kYR0lIS8/fbbSUpK4syZMxw4cIDs7GxeeOEFW/nevXuzcuVKTpw44dN98CQIMwZLGIj7ge9ExJUyWL6IyF0iskdEfhGR0R7y9RORHBG5oTDtGAwlndw4LbVq1aJnz55s2bLFtm306NE8+OCDPPXUU1SrVo3q1aszadIkbr31VptwzKeffsqRI0f46quvuP766wFLULmxY8c6aMvm0qhRIw4ePEhCQgKhoaFkZ2dz9OhR+vbtS40aNWjSpAn//ve/3do7Y8YM6tevT82aNXnllVc87luDBg0IDQ0FLDfogIAA9u3b5zJvdnY2K1ascAh9/f3339OxY0fCw8OpU6cOw4YN48qVK7btAQEBvPfeezRp0sQWbnnPnj306NGDGjVq0KxZM5sYDMCiRYu44YYbCAsLIy4ujokTJ+JrZs+eTZ8+fbjtttsICgpi0qRJfPnll5w/fz7fsp988gn9+/e3hfWuU6eOTdsgJyeHChUqOBy/ypUrc+ONN7oNf11o3IUJBXYC1ay/awLfu8vroY4AYB8Qh+Xr4S1AUxf5goFVQBJwg5u63IZBzY+SEg46L+U9HLSvuZZzpKixD8ucnJysrVq10hEjRqiq6oULF7RChQqamJjoVG7atGkaExOjqqoDBw7UBx98sMDt5oZ7Vv1NezgrK0u3bNmiNWvW1JUrV6qqczjo4OBgXbt2rWZlZekzzzzjMRy0qurs2bM1NDRURURr1aplC3Wdl9y67fnxxx9148aNmpOTo4cOHdLmzZvrm2++adsuItqjRw9NT0/XS5cu6fnz57VevXr66aefak5Ojm1fdu/eraqqq1at0h07dqiq6vbt27V27dr69ddfu7Tn8OHDWr16dQ0PD9fq1as7/A4PD3cIW21P3759dfLkyQ5pwcHB+tNPP7k9Rqqq58+f15CQECd937Vr12pYWJiKiAYHB+uyZcsctg8fPlxHjhzptl535z8ewkF7mgV0WVXPW++8J0XEm/GCvNwM7FXVQwAi8hnQF9iTJ98k4O9YxhgMhiLBV6pwhRUXuueeewCLlGD37t0LJQnZvn3BZ19rHknIJUuWOElCxufZJ3tJSIBJkybxzjvveGxn0KBBDBo0iP379zN9+nSioqJc5nMnCZmLvSTkcLuBslxJSLAI5eRKQgIOkpDjxo2jc+fOtnL2kpCuxkFyJSELii8lIcESEvzMmTMcPXqUjz76iNjYWIftISEhPg8d48kBXCciX1p/C9DQbh1VvdeL+usAyXbrR7A4BRsi0g6oq6qLRaRcOID169ezadMmKleGEiBUVibIvcl5zOPng10eJCEBGjZsSPPmzXn88cedVL3ASEK6koS0Jzo6mp49ezJw4ECH/02xSkIC/fKse3b/hUAsUw5eB4baJ7vLn/vEBBAfH+/01FJaWLRoEX/729948EF47DF/W1P6ycjIcCv4XZLIdVJlVRLSnuzsbA4cOOBym5GEzF8S0tXx2717t1cD54mJiSR629t1927IFwtwK7DEbv05YLTdeihwAsvU0oPARSy9BKdxAMrQGMBTTz2lgD75ZMl9Z11ayM7O1tmzZ+uWLVtKzRiAatmThPz3v/+tJ06csJVt0aKFPvvss26Ph5GEdGTWrFl6+PBhVVX99ddftUuXLtq/f3/b9kuXLmlERIQePXrUbd3uzn88jAEUtQOowG+DwJWwDAI385B/JdDOzTa3O54fJc0BDBkyRAEdPbrk3rBKE7k6qSXZATRo0MDp5vnEE084XOTr1q3T+Ph4DQ4O1rCwME1ISNBdu3Y5lMnIyNARI0ZovXr1NCQkRBs1aqQjR47UtLQ0r9pNSUnRhIQEjYiI0EaNGumHH35o22bvAFRVp0+frrGxsRoZGamvvPKKy33I5aGHHtKoqCgNDg7WBg0a6OjRo/Xy5ctuj8fChQu1V69etvXVq1dr06ZNNSQkRDt37qzjx493cAABAQEODkBV9ZdfftG7775ba9asqZGRkdq9e3fdunWrqqrOmzdP4+LiNDQ0VHv37q3Dhg3zuQNQVZ0zZ47GxsZqcHCw/v73v9f09HTbtl69eumrr77qkL9Zs2Y6bdo0p3rGjh2rdevW1eDgYK1Xr54+9thjDv/TuXPnar9+/TzaUuIcgKVt7gJ+BvYCz1nTJgIJLvKucPX0r2XMAfTu3VsBnTSp5N6wSiMl2QEYnLn99tt1y5Yt/jajVHDrrbfqzp07PeYpjAPwJhYQACJSWVUve5s/F1VdAlyfJ228m7zdClp/aSQ3FHRwsJ8NKWWcP3+eoKAgEz6jjGAkIb1n/fr1RVJvvlM7ReRmEdmO5QkeEWkjIm8XiTXlhD/84Q88+eST1K7tb0tKB6rKzp07+e9//1uo6XoGg8E1YukheMggsgEYAHylqu2saTtUtWUx2Gdvh+ZnqzsSJZF4jfetQT4gMVGIjy85UUpLIrkzfK5evUp8fLzHaXAiQmHPEYOhtOPu/Lemu+w2e/NxV4BaP+Syo2QE4yjhZGdnM2LECKP7Wwhyn/rnz59PbGwsffr0MeppBoOP8WYMIFlEbgZURCoAw4BfitasomHz5s289tprZGVlOaQHBgby+eefO+XPyspi4MCBTune5k9NTWXjxo0kJSWxYcMG8+66AOTk5HDq1Cn69u1rbvwGQxHhjQN4HHgLiAWOA8usaaWOtWvXkpycTFJSkkN65cqVXebPycmxhWstbP7KlSvz1ltvmZt/AalQoYLTR1EGg8G3eCMKfwJwfgwuhQwbNowePXqwc+dOh/SAANdvwgIDA11+yl6Q/DfeeCNxcXGFtNhgMBiKDm8GgT8CnDKp6p+Lyig3dphB4DKIqrJnzx6uu+46tz0rbzGDwIbyTFENAi8DlluXdUAtoMDfAxgMeckVxvj555+dxmUMvqGkSkJmZWXRokULjh8/7pP6yjJZWVk0a9asQLGYvCVfB6Cqn9stnwL3AjfmV85gcIerGT75RVAs7ZR1SchcsrOzadasmVMo47x8+OGHdOnSxW3I6NLC1KlTiY6Opnr16jzyyCNkZ2e7zTt37lyaN29OWFgYLVu25Ouvv/aqrkqVKvHwww87KMj5isLE+G8AlO7/msFvXL16lYULF7J371769u1LmzZt3I6plCXKuiRkLpMnT/bqpl4WJCGXLl3K5MmTWblyJYcOHWL//v2MH+8yyAGpqakMGTKEN954g7NnzzJ58mRbSHBv6ho0aBCffvqpRwdTKNzFiNDfYvCkA2nW5QyWyJ335VfO1ws+iAW0aNEiXbp0qccgVcVJeVUES05O1qtXr/q83ms5R4qavFE5R40apQkJCbb1XKWuvPTq1UuHDh2qqqofffSR1q5dWy9cuOBVmw0bNtQKFSpo1apVNSQkRLOysjQ1NVX79OmjERER2rhxY/3oo49s+V0Fg4uLi9PIyEh9+eWXnfYhLwcOHNDmzZvrkiVLtF69em7zHT58WIOCghzOgYULF2q7du00NDRUY2NjdcKECbZtv/76q4qIfvzxxxobG2uLprl+/Xrt2LGjVq9eXdu2beugqDZt2jRt1qyZhoSEaMOGDfWDDz7w6pgVhMGDB+vYsWNt68uXL9fatWu7zLtx40aNiopySKtZs6Zu2LDB67qaNGnipCJmj7vzn8IGg8MSmz8WS1TPClgHjf2x+MIBhIaGKuAQsc+flFcHUFSUFgdQViUhExIS9Ouvv9bExESPDmDhwoXasmVLhzRPEo65DmDo0KF64cIFvXTpkqakpGiNGjV0yZIlqqq6bNkyrVGjhp46dUpVLQ97Bw8eVFVLpNGgoCDdvHmzS3vWrl3rURJy3bp1Lsu1adNG586da1s/deqUBgQEuIzMevXqVY2Pj9f//e9/evXqVZ0/f77Wq1fP5sy9qatPnz769ttvuz2uhXEAHqeBqqqKyCIt5rAPRYGq2sSac8UwDOWLREn0ST2FnVFWliUh58+fT05ODn369MlXnMeVJGR+Eo4iwsSJE21CODNnzuTuu++mZ8+eAHTv3p327duzaNEihgwZQq9evWz1derUiR49erBmzRratm3rZM9tt93mE0nIsLAwVJXMzEzCw8Md8gYEBDBkyBAGDx7MpUuXqFy5Ml988YVtf7ypKyQkhDNnzhTYTk948yHYFhFpp6qbfdpyMZOVlcXVq1cJDAwkMDDQ3+aUeTIyMhQywpwAACAASURBVFizZg233HKLTdLQ3/h7KnBZlYS8cOECo0ePZvHixUD+8pyuJCE3bdrEc88951HCMa8k5Ny5c1mwYIGtzStXrtCtmyWg8OLFi3nppZf45ZdfyMnJ4eLFi7Ru3dqjXQUlryRkRkYGIuJyQsOyZcsYNWoUq1evpl27dvzwww/06dOHJUuW0Lp1a6/qKgpJSLejbyKS6xzaAd+LyM8i8pOIbBaRn3xqRTGQO+PCPP0XLWo3w6du3bo2zVaDa0lIwEESMi95JSGXLl3KxYsXC9W+vSRkLp4kIZOTf5Pz9iQJuXfvXg4dOkSnTp2Ijo6mX79+pKamEhMTw+HDh53yt27dmoMHD5KTk2NLGzx4MPfccw8pKSmcOXOGRx991MmR5JWEfOCBB0hLSyMtLY309HQyMzMZNWoUWVlZ9O/fn1GjRnHy5EnS09Pp1auXW8e0du1aQkJCCA0NdVhy09atW+eyXF5JyC1bthAVFeX09A+wdetWunTpQrt27QBo3749t9xyC8uWLfO6rt27d9OmTRuXthQad++GgJ+sfxu6WtyVK6qFaxwDSE5OVkDr1KlT6Hp8TVkbAzh79qz+73//0/nz5/tlnOVazpGipixLQl69elWPHz9uW7788kutU6eOnjhxQnNyclza1aZNG12/fr1t3ZOEY+4YgP2gcXJyskZHR+vSpUv16tWrevHiRU1MTNSUlBTNzMzUihUr2gZMFy1apEFBQTpu3DhP/6ICs2TJEo2OjtZdu3Zpenq6duvWTceMGeMy76pVq7RmzZo2AZyffvpJa9SoocuWLfOqrpSUFI2MjLSp37nC3flPYQaBgc3utvljuVYHkJqaqr1799YHHnig0PX4mrLkAK5evaqff/65btmypUhm+HhDSXYAZV0S0p78BoFVVd977z19/PHHbeueJBx//fVXDQgIcDqvNm3apF26dNGIiAitVauWJiQkaHJysq3+qKgoDQ8P1wceeEAHDRrkcwegqjp16lSNiorSsLAwffjhhx1u0C1atNDZs2fb1t99911t1KiRhoaGasOGDXXq1Kle1zV58mQdOXKkR1sK4wDchoIQkSPA6x56Dm63FQUmFETJ58qVK1Ss6LXInM8xoSBKD1lZWdxwww0sX7681H8MVtRkZWXRtm1bVq9e7XE8rTChIDxdrRWAYCxTQQ2GfPHnzd9QuqhUqRI7duzwtxmlgkqVKrFr164iqdvTFXtUVV8qklYNpZpz584RFBRULr7gNRjKMp6uYPPkb3BArTN85s2bVySBqQwGQ/HiqQfQvdisMJR47LV5jUqXwVA2cOsAVNV1eMFSyr59+9i1axdNmjShadOm/jan1KCq7Nq1ix9++IG2bdvSqlUr8+rHYCgjlJsr+auvvqJv3758+OGH/jal1HHmzJlyFbnTYCgvlJtpGyYOUOEQEVs8GIPBULYoN49zxgEYDAaDI8YBGIDf3vXnVakylG6MJGTZoH///ixdutTn9ZYbB2CCwbknV5t37969XLlyxd/mlEnKuiTkxIkTqVSpkkMQtVwlM1cYSUhnScgXXniBunXrEh4eTrdu3Rw+/ho9ejRjx471uf3lxgG0bt2a3/3ud9SvX9/fppQYVJUdO3Ywf/584uLi6N27N6Ghof42q0xSHiQhBw4cSEZGBpmZmWRkZHi81owkpKMk5Ny5c/nkk09Yt24daWlp3HrrrQ7H56abbiIzM5OffvJxIGZ3QYJK2oIPFMFKGv4MBnf16lX95ptv9KuvvioxCmnXyrWcI0VNWZeEzFvWE0YS0oK9JORrr72mAwYMsG3buXOnVq1a1SH///3f/+lLL73k1h535z8egsGVmx6AwZGAgABuuOEGevfubT7qKmaOHDnC4sWLady4MQAXL14kKSmJ/v37O+W97777+O677wBYvnw5d911l01FKj/27dtHvXr1bD2PwMBABgwYQGxsLMeOHeOLL75gzJgxJCYmOpXdtWsXTzzxBLNmzSI1NZXTp0+TkpLisb0FCxYQGRlJq1ateP/9993m2759O9ddd53DlOLg4GBmzJjB2bNnWbhwIe+//77TmMXq1avZs2cPS5cuJTU1lYSEBF588UXS09OZMmUK/fr1s32hHhUVxaJFi8jIyGDatGmMGDGCLVu2uLRn3bp1hIeHExERQXh4uMPviIgIkpKSXJbbuXOnQ3z+Nm3acOLECZfqYu3bt6dZs2YsWLCAnJwcvvrqK6pUqWITqRk4cCD79+9n7969ZGdn88knnziomgE0a9bMQTPAF5SbaaAGZ1xJEJZlEhN9E92ksBFcy7Ik5IABA3j00UeJiopiw4YN9OvXj/DwcAYMGOCU10hCOktCRkdHc9ttt3H99ddTsWJF6tWr5zTg7i9JyGtCRO4C3sAy3vCxqr6WZ/sI4BEgGzgJ/ElVk50qMhQaVS304F5Zwt+ht8uqJCTg8HV9hw4d+Mtf/sJ///tflw7ASEI6S0JOnDiRH374gZSUFKKiopgxYwZdu3bl/9s78/Cqqqvh/1YghAQCJNwwVIYwBSLKFOENKEiFj1HEFhCRSYoItNgXjV+l+lVN+7WgAvpaVDRqGFrBymBBCcWHKVRlktFARRSZApSQhDAmwl3vH/fk9t7kJmS4N8O9+/c858k9+6x99trn3uy1p7PWoUOHqF27NlDBISG9gYgEAQuAgUBHYIyIFPTDsAeIU9UuwErgFV/qFGjk7/A5e/ZsZasS8OT3xP0tJKQniovNYEJCFg4JuX//fkaPHk3Tpk0JCgpi4sSJZGVlue0EqtCQkN44gHggxeV8FvBMMfJdgG1FXCty8eNWbGazrl69Wj/55JNiQ6pVNL5cBLbb7Xrw4EFdtGiR7t+/v9KidFUk5fmN+Bp/Dgmpqvr3v//duZlgx44detttt+nSpUuLfB4mJKR7SMjExETt3bu3njt3Tu12uy5ZskTr1q2rFy9edN4jJiZGd+3aVaQ+Rf3+KUtISG8cwAjgHZfzccDrxcj/GXi2iGtFVvxWbGazhoaGKqCXL18u8328ja8MQH5sXn/a4VMSqrIB8PeQkGPGjNGGDRtqeHi4xsbG6oIFC4p9HiYkpHtIyOvXr+uMGTO0adOmWr9+fY2Li9MNGza41TUuLq5YXcpiAIoMCekNRGQEMFBVH7fOxwE9VPXXHmTHAb8E7lXVQm9TiIi67rHt27dvoYWrotgkm+hnebe+efNmlXFo5ouQkKrKqlWraNeuHXfccUeVqWtFYEJCVh9MSMjSMXLkSB577DEGDRpUpEz+73/Lli1uO7sSExPRIkJC+toAxAMvquog63wWDmtUcCG4P/A/QB9V9TjRWJ6YwCmSwhCGEBYW5jb/Wdn4Kiaw3W4PqIY/H2MADIFMWWIC+7qV2AW0FZGWIlILeBhw29wrIl2BhcADRTX+5eUajkWzQHEDEYiNv8FgKD0+bSlU9SYwA9gApAHLVfWwiCSKyP2W2MtAHeAjEdkrIh97W4/rXAf8zwDk5ORUmdfiDQZD9cPn7wGo6nqgfYG0F1w+/x9f6xBMMCNGjHDup67uqDpi83711VcMHDiQJk2aVLZKBoOhGhIQbwJHEcWKFSsqWw2vkB+b1263m9i8BoOhXASEAfAHXHv9Xbt2DbgdPgaDwfsYA1CNuHr1qun1GwwGr2EMQDVBROjRo0dlq2EwGPwIM4dgMFQzjh8/TlBQkJsvHYOhLASEATjBCVauXOnmWKmqkj/X7+pl0FD98WY8XSh76EaDwZWAMADb2MbIkSNZsmRJZatSLPmeO48ePWreaDUYDD4nIAxA/otgrn7QqxKqhWPzugaaMPgnixcvpnfv3m5pQUFBfP/99wBcv36dhIQEoqOjiYiIoE+fPuTm5ha6z8qVK2ndunW1GOEaqhYBsQhcld8EVlXWr19PXl6e2eETgBScynE9T0hI4PDhw2zfvp3GjRuzY8eOQlt/k5OTmT17Nhs3bqRVq1YVorPBfwgIA1CVfQGJCHFxcdhsNrOv38fs3r2bPXv2FErv1q2bx1CLnuSLkvUW+VN/qkpycjI7d+50vukdHx/vJvfqq6+SnJzM1q1bAy68p8E7BIQBcI4Apk+H6dMrWRsXNjv+NGrUqHL1CBDuuuuuUjXepZX3JhkZGeTm5tK6desiZebOncvzzz9vGn9DmQmILufGgRt58MEHabVtGzii4FTKoXa7e5ohoKlTpw5Xr151nruG7bTZbNSuXZvvvvvOY14RYcOGDfzhD39g1apVPtfV4J8EhAGgJ6xevZp77rmn0lTI3+HjGmfVEFjk5eWRm5vrPDp37kxaWhoHDhwgNzeXxMRE5xqAiDBp0iSeeuopzpw5g91uZ/v27fz4oyNWkqrSsWNH1q9fz4wZM5zB0Q2G0hAYBqASKbjDx1MAbkNgMHToUMLCwggNDSUsLIylS5fy/PPP069fP2JiYgrtCJo7dy533nkn3bt3p2HDhsyaNcv58le+oejUqRNr167l8ccf5x//+EeF18lQvfFpRDBvUp6IYJWFq+fOe++9t9AOH19FBAtUTEQwQyBTlohgAbEIXBmoKps3b6ZVq1bGc6fBYKiSGAPgI0SEYcOGmYbfYDBUWQKidfrb3/7Gxx9/zI0bNyq0XNP4GwyGqkxAjADGjh3LjRs3PL5G7w1ycnIIDQ0lODjYJ/c3GAwGX+D3XdS8vDxu3LhBzZo1qVWrllfv7brD59///rdX720wGAy+xu9HAFeuXAG87wbCxOY1GAzVHWMASomJzWswGPwFYwBKiYjw448/ml6/wWCo9vh91zU0NJSHH36YIUOGeO2eXbt2NY2/wVBGMjIyiI2N9dmmDH8iLy+P2NhYLly44JP7+70BaNGiBcuWLeO1116rbFUMAUx0dDS1a9cmMzPTLb1r164EBQVx4sSJCtVn69at1KhRg3r16lG/fn1iY2NZtGhRIblXXnmFmJgY6tSpQ3R0NM8++yx5eXluMjt37mTo0KFERERgs9mIj4/3eK985syZw6RJkwgJCfFyrSqWZ555BpvNRlRUFLNmzSpSbvbs2YSHh1OvXj3q1atHWFgYNWvWdP4WsrKyGD16NDabjUaNGjF+/HguX74MQK1atZg8eTKzZ8/2TSVUtVocDlUrDrvdrgcPHtSMjAyflbF5c8XWyd+p6N9IaYiOjtYOHTroggULnGkHDx7U9u3ba1BQkB4/frxC9dmyZYs2b97ceb5u3TqtWbOmHjlyxJk2Y8YMjYmJ0R07dujNmzf10KFD2qNHDx0+fLhT5osvvtC6devqK6+8ohcuXFBV1T179ujo0aM9lpubm6s2m01Pnz5dJr1v3LhRpnzeZuHChdqhQwdNT0/X9PR0vf322/Xtt98uUd4XX3xR+/Xr5zyfPn26Dhw4UC9fvqw5OTnav39/TUhIcF4/deqU2mw2zcvLK/a+Rf3+rXTP7WpRF6raUZH/3BcvXtQ1a9boxx9/rNnZ2T4rxxgA71LVDcAf//hH7d69uzPt6aef1j/96U9uBiA3N1cTEhK0RYsW2qRJE50+fbpev35dVVWzsrL0/vvv16ioKI2MjNT7779fT5065bxf37599Xe/+53efffdGh4ergMHDnQ2ygUpaABUVRs1aqQrVqxQVdUjR45ojRo1dPfu3W4yJ0+e1JCQEN28ebOqqt5zzz36xBNPlPg5pKamart27dzSkpOTNTY2VsPDw7VNmzZuDemWLVu0WbNm+tJLL2mTJk10woQJqqq6du1a7dKlizZo0EDvvvtuPXDggDPPnDlztE2bNhoeHq4dO3bU1atXl1i/ktKrVy9NSkpynr/33nvas2fPEuVt3bq1Ll261Hk+ePBgfeutt5znb7zxhg4aNMgtT0xMjKamphZ737IYAL+fAioNqiY2rz8jIh6P0siXh/j4eC5dusQ333yD3W7nww8/ZNy4cW4OvJ555hmOHj3KgQMHOHr0KKdPn+b3v/89AHa7nV/84hecPHmSEydOEBYWxowZM9zKWLZsGYsXL+b8+fPk5uYyd+7cW+qlqqxZs4YLFy7Qtm1bADZt2kTz5s2Ji4tzk23WrBnx8fF89tlnXLt2jS+//JIRI0aU+BkcPHiQ9u3bu6U1btyYdevWkZOTQ3JyMk8++ST79u1zXj979izZ2dmcOHGCd955h7179zJ58mSSkpLIzMxk6tSpPPDAA05X2W3btuXzzz8nJyeHF154gXHjxnHu3DmP+ixbtoyIiAgiIyOJiIhw+xwZGcmpU6c85ktLS6Nz587O83zX3rciNTWV8+fP8/Of/9yZ9qtf/Yq1a9eSnZ1NVlYWK1euLLRm2aFDB/bv33/L+5cWYwBc2LBhA9999x3Dhw+nU6dOZnunweuMHz+exYsX89lnnxEbG8tPfvITt+tJSUm8+uqr1K9fnzp16jBr1iyWLVsGQGRkJD/72c8ICQmhTp06/Pa3vyU1NdUt/6RJk2jTpg0hISE89NBDbg1pQU6fPk1kZCShoaGMGDGC+fPnOxu1jIyMIiONNW3alIyMDLKysrDb7aWKSJadnU14eLhb2uDBg4mOjgagd+/eDBgwgG3btjmv16hRg8TERIKDgwkJCSEpKYlp06Zx1113ISKMHz+ekJAQtm/fDsCIESNo3LgxAKNGjaJdu3bs3LnToz5jxowhKyuLzMxMsrKy3D5nZmbSrFkzj/kuX77s1jmsX7++c96+OJYsWcLIkSMJCwtzpnXr1o28vDwaNmxIVFQUNWvWZHqByIXh4eFkZ2ff8v6lxe9buIMHD7Js2bISWefu3bszbNgws8PHTylqGFwa+fIybtw4PvjgAxYtWsSECRPcrp0/f56rV68SFxdHZGQkkZGRDB482LkD5Nq1a0ydOpXo6GgaNGjAvffeS3Z2tpte+fGDAcLCwoptlG677TYyMzO5dOkSv/71r9m0aZPzms1m48yZMx7znTlzBpvNRkREBEFBQUXKeSIiIoJLly65paWkpNCzZ08aNmxIREQEKSkpZGRkOK9HRUW5uVk5fvw48+bNcz6jiIgITp06RXp6OuBoZLt27ers0aelpbndzxvUrVuXnJwc53lOTg5169YtNs+1a9f46KOPePTRR93SR40aRfv27bly5Qo5OTm0bt2asWPHuslcunTJJ+2S3xuA1atX88gjj7B8+fJbykZGRppev8GntGjRglatWpGSkuI2DQCORjcsLIy0tDQyMzPJzMwkOzubixcvAjBv3jy+/fZbdu3aRXZ2trP3X17DFBwczJw5czhw4ABr1qwB4L777uPkyZPs3r3bTfbkyZNs376d/v37ExoaSs+ePVm5cmWJy+rUqRNHjhxxnufl5TFy5Eh+85vfcP78ebKyshg8eLBbnQpOvTVv3pznnnvO+YyysrK4fPkyo0eP5sSJEzz++OO8+eabzh59x44di3xGH3zwgdsOnfwjP62oKaCOHTu6Tcns27ePjh07Flv3VatW0bBhQ/r06eOWvn//fqZOnUrt2rUJCwtj2rRppKSkuMkcPnzYbcrJW/i8tRORQSLyLxE5IiLPeLheS0SWi8i3IvKliLTwZvmeXgTzVm/OYCgL77//Pps2bSI0NNQtXUSYMmUKM2fO5Pz584BjmmbDhg2AoxcYGhpKvXr1yMzM5MUXX/SaTsHBwSQkJJCYmAhAu3btmDp1KmPHjmXHjh3Y7XbS0tIYOXIkAwYM4Kc//SkAL7/8MosWLWLevHnObY379+9nzJgxHsvp0aMH2dnZzlFDXl4eeXl52Gw2goKCSElJcda3KKZMmcLChQud0zpXrlxh3bp1XLlyhStXrhAUFITNZsNut5OcnMzXX39d5L0eeeQRLl26RE5OjtuRn1bUFNCECROYP38+6enppKenM3/+fCZNmlSs3kuWLCk06st/Ju+++y7Xr1/n2rVrvP3223Tq1Ml5PT09naysLOLj44u9f1nwqQEQkSBgATAQ6AiMEZEOBcQmA5mq2g54DXjZmzrkG4D8Obf82LzHjh3zZjEGQ7G49mJbtWpFt27dPF576aWXaNu2LfHx8TRo0IABAwY4e8wzZ87k6tWr2Gw2evXqVWihsLyL1PkLzJ9++ikAb7zxBo899hjjxo0jPDycIUOGcN9997FixQpnnp49e7Jp0yY2btxImzZtsNlsTJs2jaFDh3osIzg4mEcffZSlS5cCjqmU119/nVGjRhEZGcny5csZPnx4sXrGxcWRlJTEjBkziIyMJCYmhsWLFwMQGxtLQkIC8fHxNGnShLS0NJ/EAp86dSrDhg3jzjvvpFOnTgwbNowpU6Y4r4eHh/P55587z9PT09m8ebNHA/D+++9z7NgxmjVrRvPmzfnhhx+c9QH461//ysSJE33ibdinISFFJB54QVUHW+ezcGxJeslFZr0ls0NEagBnVTXKw720LLpOmjSJRYsWkZSURHx8fJXy4WNCQnoXExKyepCRkUGfPn3Yu3dvtX8ZzNfk5eXRpUsXUlNTsdlsxcpWxZCQtwEnXc5PAT2KklHVmyKSLSKRqpqJF8gfAXzzzTdERUUZHz4GQyVjs9k4dOhQZatRLahVq5ZPn1VVdAZXvnFsAXr16sXZs2fp0qWLCdFoMBgMLvjaAJwGXBd1m1lprpwCmgPp1hRQvaJ6/66LXn379qVv3763VGDmzJl07vwUIttITR1XKuV9Tc2aEZWtgsFg8DO2bNnCli1bSiTr6zWAGsA3QD/gDLATGKOqh11kfgncoaq/FJGHgQdV9WEP9yrTGoAhcDBrAIZApsqtAVhz+jOADTh2HL2nqodFJBHYpaqfAO8BS0XkW+ACUKjxNxgMBoP38ekIwJuYEYDhVpgRgCGQqXIjAIOhImnZsmW598IbDNWVli1bljqPGQEYDAaDH1PcCCAg9kSWdEXcnzB1DgxMnQMDX9XZGAA/xdQ5MDB1DgyMATAYDAaDVzEGwGAwGAKUarUIXNk6GAwGQ3WkqEXgamMADAaDweBdzBSQwWAwBCjGABgMBkOA4lcGoLLDT1YGJajzkyKSJiL7ROQzEWleGXp6k1vV2UVuhIjYRaRbUTLVhZLUWUQesr7rgyLyl4rW0duU4LfdXEQ2icge6/c9uDL09BYi8p6InBORA8XIvG61X/tEpEu5C82Pj1vdDxzG7CjQEggG9gEdCshMB960Po8Glle23hVQ53uB2tbnaYFQZ0uuLrAV+ALoVtl6V8D33Bb4Coc7dQBbZetdAXV+G5hqfY4FjlW23uWs8z1AF+BAEdcHA59an/8L2F7eMv1pBNAD+FZVj6vqj8ByoGBw0eFAfrDNFTjcVFdnbllnVd2qqtet0+04IrBVZ0ryPQP8AZgD5Fakcj6iJHWeAryhqjkAqppRwTp6m5LU2Q7Usz43oHCskWqFqv4TyCpGZDiwxJLdAdQXkcblKdOfDICn8JMFGzu38JNAtohEVox6PqEkdXZlMpDiU418zy3rLCJdgWaqWt3rmk9JvucYoL2I/FNEvhCRgRWmnW8oSZ0TgfEichL4BHiignSrLAo+k9OUs0MX6N5AA8Z1pIiMA+JwTAn5LeJwBzofmOiaXEnqVCQ1cUwD9cERhS9VRO7IHxH4KWOAZFV9VUTigb8AHStZp2qFP40AShN+Mj9aWZHhJ6sJJakzItIf+C0wzBpOV2duVedwHI3AFhE5BsQDf6/mC8El/W2vUVW7qv4AHAHaVYx6PqEkdZ4M/A1AVbcDtUXEVjHqVQqnsdovC4//76XBnwzALqCtiLQUkVo4IoutKSCzlv/0DEcBmypQP19wyzpb0yELgQdU9UIl6Ohtiq2zquaoaiNVba2qrXCsewxT1T2VpK83KMlv+2PgpwBWI9gO+L5CtfQuJanzcaA/gIjEAiF+sPYhFD1iXQNMALBGPNmqeq48hfnNFJAGYPjJEtb5ZaAO8JE1PXJcVR+sPK3LRwnr7JaFaj4FVJI6q+o/RGSAiKQBN4CnVbW4BcUqTQm/56eBJBF5EseC8MSi71j1EZEPgL5AQxE5AbwA1AJUVd9R1XUiMkREjgJXgEnlLtPaUmQwGAyGAMOfpoAMBoPBUAqMATAYDIYAxRgAg8FgCFCMATAYDIYAxRgAg8FgCFCMATAYDIYAxRgAQ5VBRG5arn33Wn+LdNdtvSB00AtlbrZcDu8TkW0iUuq3Z0VkquVqAxGZKCJNXK69IyIdvKznDhHpVII8/y0itctbtsF/MQbAUJW4oqrdVLWr9ffELeS99RLLGFXtgsPT4tzSZlbVt1U13//+o7g46FLVx1X1X17R8j96vkXJ9JwJhHmpbIMfYgyAoSpR6I1dq6efKiK7rSPeg8ztVq84PzBIGyt9rEv6W9ab0MWVmwrk5+1n5dsvIu+KSLCVPkdEvrbKedlKe0FEEkRkBHAX8Bcrb22r597NGiW87KLzRBF5vYx6fgn8xOVeb4rITnEEgnnBSnvCktksIhuttAGWp9DdIvKhiBjjEOAYA2CoSoS6TAGttNLOAf1V9S4crjv+7CHfNOA1Ve2GowE+ZU27jAZ6Wel2YOwtyn8AOCgiIUAyMEpVO+MISDJdHK7DH1TVO6ye+P93yauquhLYDTxijWCuu1xfCfzM5Xw0sLyMeg7C4fsnn2dVtQfQGehreQH9Mw5HYX1VtZ+INASeA/pZz/IrIOEW5Rj8HL/xBWTwC65ajaArtYAF4gh/dxPPHi6/BJ4TR7jLVap6VET6Ad2AXVaPujYOY+KJv4rINeAHHD7l2wPfq+p31vXFwC+BN4BrIvIu8CkOH/SeKNSDV9UMEflORHrgiHTVXlW/EJFflVLPEBy+nVzDAT4sIlNw/D83AW4HvsbdsVi8lf65VU4wjudmCGCMATBUdZ4EzqpqJ3G48L5WUEBVl4nIduB+4FMRmYqj4Vusqs+VoIxHVHVv/onVW/bUiN+0GvB+OLzJzqB0UeU+xNHb/xewOr+40uppaKugyQAAAWVJREFUTSUtAEaISDSOnnycquaISDIOI1IQATao6q1GF4YAwkwBGaoSnua+6wNnrM8TgBqFMom0UtVj1rTHGqATsBEYKSJRlkxEMbuKCpb7DdBSRFpb5+OBrdaceQNVXQ88ZZVTkEv8J0xhQVbjCOv3MI4Qh5RRz+eB/xKRGKusy8AlcYQHdA2MnuOiy3bgbpf1kbCy7Hgy+BfGABiqEp529bwJPCoie3GEPbziQeYha2F2L45gMEtU9TDw/4ANIrIfh1vhJh7yFipTVXNxuNpdYeW9iSOmQj3gEystFcfopCCLgIX5i8Cu91fVbOAw0EJVd1tppdbTWluYB/xfVT2AI2D6YRwRsf7pkicJWC8iGy0/+ZOAZVY5X+CY6jIEMMYdtMFgMAQoZgRgMBgMAYoxAAaDwRCgGANgMBgMAYoxAAaDwRCgGANgMBgMAYoxAAaDwRCgGANgMBgMAYoxAAaDwRCg/C9j9fF2dMDVAgAAAABJRU5ErkJggg==\n",
      "text/plain": [
       "<matplotlib.figure.Figure at 0x7f65e2259710>"
      ]
     },
     "metadata": {},
     "output_type": "display_data"
    }
   ],
   "source": [
    "# -*- coding: utf-8 -*-\n",
    "\"\"\"\n",
    "Created on Sun Apr 19 08:57:13 2015\n",
    "\n",
    "@author: shifeng\n",
    "\"\"\"\n",
    "print(__doc__)\n",
    "\n",
    "import numpy as np\n",
    "from scipy import interp\n",
    "import matplotlib.pyplot as plt\n",
    "\n",
    "from sklearn import svm, datasets\n",
    "from sklearn.metrics import roc_curve, auc\n",
    "from sklearn.cross_validation import StratifiedKFold\n",
    "\n",
    "###############################################################################\n",
    "# Data IO and generation,导入iris数据，做数据准备\n",
    "\n",
    "# import some data to play with\n",
    "iris = datasets.load_iris()\n",
    "X = iris.data\n",
    "y = iris.target\n",
    "X, y = X[y != 2], y[y != 2]#去掉了label为2，label只能二分，才可以。\n",
    "n_samples, n_features = X.shape\n",
    "\n",
    "# Add noisy features\n",
    "random_state = np.random.RandomState(0)\n",
    "X = np.c_[X, random_state.randn(n_samples, 200 * n_features)]\n",
    "\n",
    "###############################################################################\n",
    "# Classification and ROC analysis\n",
    "#分类，做ROC分析\n",
    "\n",
    "# Run classifier with cross-validation and plot ROC curves\n",
    "#使用6折交叉验证，并且画ROC曲线\n",
    "cv = StratifiedKFold(y, n_folds=6)\n",
    "classifier = svm.SVC(kernel='linear', probability=True,\n",
    "                     random_state=random_state)#注意这里，probability=True,需要，不然预测的时候会出现异常。另外rbf核效果更好些。\n",
    "\n",
    "mean_tpr = 0.0\n",
    "mean_fpr = np.linspace(0, 1, 100)\n",
    "all_tpr = []\n",
    "\n",
    "for i, (train, test) in enumerate(cv):\n",
    "\t#通过训练数据，使用svm线性核建立模型，并对测试集进行测试，求出预测得分\n",
    "    probas_ = classifier.fit(X[train], y[train]).predict_proba(X[test])\n",
    "#    print set(y[train])                     #set([0,1]) 即label有两个类别\n",
    "#    print len(X[train]),len(X[test])        #训练集有84个，测试集有16个\n",
    "#    print \"++\",probas_                      #predict_proba()函数输出的是测试集在lael各类别上的置信度，\n",
    "#    #在哪个类别上的置信度高，则分为哪类\n",
    "    # Compute ROC curve and area the curve\n",
    "    #通过roc_curve()函数，求出fpr和tpr，以及阈值\n",
    "    fpr, tpr, thresholds = roc_curve(y[test], probas_[:, 1])\n",
    "    mean_tpr += interp(mean_fpr, fpr, tpr)\t\t\t#对mean_tpr在mean_fpr处进行插值，通过scipy包调用interp()函数\n",
    "    mean_tpr[0] = 0.0 \t\t\t\t\t\t\t\t#初始处为0\n",
    "    roc_auc = auc(fpr, tpr)\n",
    "    #画图，只需要plt.plot(fpr,tpr),变量roc_auc只是记录auc的值，通过auc()函数能计算出来\n",
    "    plt.plot(fpr, tpr, lw=1, label='ROC fold %d (area = %0.2f)' % (i, roc_auc))\n",
    "\n",
    "#画对角线\n",
    "plt.plot([0, 1], [0, 1], '--', color=(0.6, 0.6, 0.6), label='Luck')\n",
    "\n",
    "mean_tpr /= len(cv) \t\t\t\t\t#在mean_fpr100个点，每个点处插值插值多次取平均\n",
    "mean_tpr[-1] = 1.0 \t\t\t\t\t\t#坐标最后一个点为（1,1）\n",
    "mean_auc = auc(mean_fpr, mean_tpr)\t\t#计算平均AUC值\n",
    "#画平均ROC曲线\n",
    "#print mean_fpr,len(mean_fpr)\n",
    "#print mean_tpr\n",
    "plt.plot(mean_fpr, mean_tpr, 'k--',\n",
    "         label='Mean ROC (area = %0.2f)' % mean_auc, lw=2)\n",
    "\n",
    "plt.xlim([-0.05, 1.05])\n",
    "plt.ylim([-0.05, 1.05])\n",
    "plt.xlabel('False Positive Rate')\n",
    "plt.ylabel('True Positive Rate')\n",
    "plt.title('Receiver operating characteristic example')\n",
    "plt.legend(loc=\"lower right\")\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": []
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.5.2"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
